Mechanistic model of the Oxygen isotope in tree ring cellulose

Isotope Dendrochronology

COURSE: DISC, LTRR, UNIVERSITY OF ARIZONA. INSTRUCTORS: P. SZEJNER & S. BELMECHERI

Instructions :

The following code and comments should further your experience of using a mechanistic model to understand the isotopic fractionation mechanisms underlying Oxygen isotope (d18O) in tree ring cellulose. This exercise will focus on the main determinants (variables and parameters) of cellulose d18O including: (1) source water isotope ratio, (2) leaf water isotopic enrichment, and (3) biochemical fractionation at the site of sucrose production and cellulose synthesis.

The exercise and the mechanistic model correspond to the process based models termed “Craig-Gordon” model [1], and the “Roden” model [2].

[1] Craig H, Gordon L (1965) Deuterium and oxygen 18 variations in the ocean and the marine atmosphere. In: Tongiorgi E (ed) Stable Isotopes in Oceanographic Studies and Paleotemperatures. Consiglio Nazionale delle ricerche Laboratorio di Geologia Nucleare. Pisa, Spoleto, Italy, pp 9–130).

[2] Roden J, Lin G, Ehleringer JR (2000) A mechanistic model for interpretation of hydrogen and oxygen isotope ratios in tree-ring cellulose. Geochim Cosmochim Acta 64:21–35. https://doi.org/10.1016/S0016-7037(99)00195-7

The stable Oxygen isotope composition of cellulose (δ18Ocell) is determined by three factors/processes:
1. The isotopic content of the source water.
2. The leaf evaporative enrichment effect.
3. The proportion of phloem sugars that exchange with source water during cellulose formation.

The signal of these three processes is integrated into the cambium during xylogenesis. The goal of this exercise is to: (1) implement a mechanistic model to account for all these factors processes and their associated fractionations, (2) to use this mechanistic model to interpret variations of cellulose d18O.

You will simulate three different isotopic ratios in tree ring cellulose for different environmental conditions. These environmental condition reflect the growing season variability of relative humidity and precipitation experienced by trees growing in the US Southwest (See Figure 1). The growing season in this region is characterized by a bimodal distribution of precipitation (Winter-snowpack and Summer-Monsoon). The relative humidity varies from dry spring to humid summer. The two periods are separated by a hyper arid period with low precipitation and extremely low relative humidity. Accordingly, the variables in the model will vary as follows:

• The relative humidity from 50% (spring) to 30% (pre-summer drought) to 60% (monsoon).
  
• The mean air temperature will vary from 25, to 33, to 27 Celsius.
  
• The source water will remain constant at -8 permill.

At the end of the exercise, you will compare the simulated 18O cellulose for the 3 seasonal variations of RH (and temperature) with measurements of cellulose d18O in tree rings from a site located in the US Southwest.

Figure 1: Seasonal relative humidity in HEL-Catalina Mountains (Arizona) (Belmecheri et al. 2018 PCE). Notice the large variations of relative humidity.

• You will need to download the R project from https://github.com/paulszejner/DISCO-Notebook2.git and set your own working directory. Once you open the Notebook2.Rproj in R studio you can open the Notebook2.Rmd file to proceed with this exercise.

PART 1: Source Water

The d18O of source water noted as Source_18O is the isotopic ratio of the water transported through the xylem up to the leaves. This water is derived from direct precipitation and soil water taken up by the roots. This process does not alter the original d18O. Thus, source water is the starting point of d18O that will be modified by the isotopic effects (fractionations) during all the steps leading to its fixation in tree ring cellulose. Having a robust estimation or measurement of source water d18O is therefore paramount.

Here, the source water is set up at a fixed value:

Source_18O <- -8 # This is the isotopic value of the source water and  will remain constant at -8 permill during this  exercise

Then, it is converted to an isotopic ratio with regard to the standard VSMOW:

rsmow <- 0.0020052 # Oxygen Isotope ratio of the STANDARD   VSMOW
Rsource= rsmow*(1+Source_18O/1000) 

PART 2: leaf water isotopic enrichment

During plant transpiration, the heavier H218O evaporates and diffuses more slowly through the stomata than does H216O, leaving leaf water enriched in 18O. At steady state and considering the leaf as a single and isotopically well-mixed water pool, the leaf water isotopic enrichment can be estimated using the Craig Gordon model as follows:

• De is the isotopic enrichment of evaporative site water (leaf) above the source.
  
• Epsilon_plus is the equilibrium fractionation between liquid water and vapor.
  
• Epsilon_k is the kinectic fractionation for combined diffusion through stomata and the boundary layer.
  
• Dv is the isotopic enrichment of atmospheric vapor compared to source water.
  
• wa/wi is the ratio of the water vapor mole fraction in the air relative to that in the intercellular air spaces. This is also relative humidity assuming leaf and air temperature are equal.

Climate Variables/Inputs

# The 3 temperature  and Relative humidity values need to be modified to reflect the seasonality described in the introduction

temperature <- c(40, 40, 40) # temperature here is in Celsius 
Relative_humidity <- c(100, 100, 100) # RH is in %

Constants

These constants will be used to estimate/calculate the various model parameters, specifically the water vapor mole fraction in the air and in the intercellular air spaces.

br <-1 # boundary layer conductance. (mol m-2 s-1), considered reasonable value between 0.02 to 3 in Roden paper boundary layer resistance.
ak <- 1.0285    # alpha_kappa, Kinetic Fractionation 
akb <- 1.0189   # alpha_kappa_beta, Kinetic Fractionation boundary layer
Stomatal_conductance <- 0.5 # stomatal conductance to water vapor
press_bar <- 77 # barometric    pressure  (mmbar)  

Model Parameters/Variables

Water Vapor (Dv)

The isotopic value of the water vapor is estimated based on source water d18O and leaf temperature, assuming equilibrium between atmospheric vapor and source water. Water vapor can be estimates as follows:

water_vapor = Source water + the equilibrium fractionation between liquid water and vapor.

The equilibrium fractionation varies as a function of temperature.

# the temperature in the equation below refers to the leaf temperature in degrees Celsius and is converted to Kelvin (273.15)
water_vapor <-Source_18O +(-1*(1.137*(1000000/((temperature+273.15)^2))-0.4156*((1000/(temperature+273.15))-2.0667))) 

Then, it is converted to an isotopic ratio with regard to the standard VSMOW:

Rwater_vapor= rsmow*(1+water_vapor/1000)
#note!  I need help on figuring out this! What is needed?

Alpha *

The equilibrium fractionation at leaf temperature is temperature dependent (see Majoube, 1971, J. Chim. Phys. 58:1423-1436)

alpha_star =exp((1.137*1000/(273.15+temperature)^2)-(0.4156/(273.15+temperature))-0.0020667)

Water vapor mole fraction in the intercellular air spaces

Here we assume that the vapor is saturated inside the leaf. We apply this formula to calculate the vapor pressure when saturated (Rh=100%) at the temperature we best estimate the plant is and call it "ei".

ei <- (101325*exp((((-0.1299*(1-(373.16/(273.15 + temperature)))-0.6445)*(1-(373.16/(273.15 + temperature)))-1.976)*(1-(373.16/(273.15 + temperature)))+13.3185)*(1-(373.16/(273.15 + temperature)))))/1000

Then leaf water vapor is converted to mole fraction as follows:

wi=ei/press_bar

Water vapor mole fraction in the air

ea <- (Relative_humidity/100)*ei

Then water vapor in the air is converted to mole fraction as follows:

wa= ea/press_bar

Total Conductance

g = 1/(1/Stomatal_conductance + 1/br) # or 1/(rs+rb)  in mol

Leaf transpiration

Transpiration will depend on water vapor mole fraction in the air, in the leaf and on total conductance as follows:

E <- (wi-wa)*g #(mol m-2 s-1) # barbour 2004

Leaf surface water vapor

ws <- ((Stomatal_conductance*wi)-E*(1-wi/2))/(Stomatal_conductance - E/2) #mole fraction Roden (2000), Ball (1987).

Vapor pressure at leaf surface

esl <- ws*press_bar #(kPa) higher pressure  than the atmosphere

Estimating Leaf water isotopic enrichment

Lets recall the estimated variables necessary to calculate leaf water

Rsource Source water
Rwater_vapor isotopic value of the water vapor
ei Water vapor pressure inside the leaf
ea Water vapor pressure of the air
esl Water vapor pressure at leaf surface

The values can be rounded as follows:

round (Rsource, digits = 4)
round(Rwater_vapor, digits = 4)
round(esl, digits = 2)
round(ea, digits=2)
round(ei, digits=2)

The leaf water enrichment can be estimated as follows:

Rleaf <- alpha_star*((ak*Rsource *
                              ((ei-esl)/ei))+
                             (akb*Rsource *(esl-ea)/ei)+
                             (Rwater_vapor*ea/ei))

The leaf water enrichment is estimated in terms of isotopic ratio (R). It should be converted to permill values as follows:

rsmow <- 0.0020052
leaf_water_18O <- ((Rleaf/rsmow)-1)*1000
round(leaf_water_18O, digits = 2)

PART 3: isotopic signal in the Cellulose

Isotopic fractionation during sugar synthesis and mixing during cellulose synthesis.

Now, after we know the water where the sugars are made and transported in the phloem, we can use the leaf water and the Source_18O to estimate the result on the cellulose, take into consideration that during cellulose synthesis, the source water will have some influence on the result, and we usually use a 42 percent exchange.

# Now   after we know the  water where the sugars are made and trasported in the floem  we can use the leaf water and the Source_18O to estimate the end result on the cellulose, take in consideration that during  cellulose sinthesisi  the source water will have some influence on the end result. and we usualy use a 42 percent exchange
# and some constants

epsilon_s <-  27 #Autotrophic fractionation Sucrose 
epsilon_c <- 27 # Autotrophic fractionation Cellulose ## 0.0084*Tmin.ts^2-0.51*Tmin.ts+ 33.172 
F_0 <- 0.42 # proportion of source water used during cellulose sinthesys

cellulose_d18O= F_0*(Source_18O + epsilon_c) + (1 - F_0) * (leaf_water_18O+ epsilon_s) 

PART 4: Comparison between modeled data and observations

Mechanistic modeling of tree‐ring cellulose δ18O and leaf water enrichment ΔL, and their sensitivities to relative humidity. (a) The Craig–Gordon model (CG) outputs are plotted from May to December, the period corresponding to cell wall thickening for all tree‐ring subdivisions. The yellow color lines correspond to the CG simulations of δ18O cellulose within the range of the proportion of exchangeable oxygen in cellulose (Pex) variability from 0.1 (thin line) to 0.7 (thick line). The gray horizontal bars correspond to the cell wall thickening (or cell‐maturation) period for each subdivision and are plotted on the mean of the δ18O measurements. The cell wall thickening phase was derived from microscopic phenological observations from Morino and Hughes (University of Arizona, unpublished observations, Table S3). The climate drivers used in the CG model to make δ18Ocel predictions were derived from PRISM data; (b) seasonal cycle of the CG mean leaf water enrichment (ΔL). The shaded gray color corresponds to the standard deviation of modelled ΔL. The dashed blue line in (a) and (b) indicates the average date of the monsoon onset; (c) modelled ΔL plotted versus relative humidity (RH) in % calculated at bimonthly resolution. The ΔL simulations were run using mean, minimum, and maximum RH. The range of RH for the hyper‐arid period of May and June are plotted in orange circles. The dark blue circles represent RH range observed for the monsoon season (July–August). The gray circles correspond to the range of RH observed from September to April (Belmecheri et al. 2018 PCE)

your estimations for the callulose oxygen Isotopes are:

round(cellulose_d18O, digits = 2)

your estimations for the leaf water oxygen Isotopes are:

round(leaf_water_18O, digits = 2)

Comparison between modeled data and observations,

Now with these values that you comared write down your observations differences and similarities between the modeled data in this excersise and the measurements shown in the figure 2 above.

How do you think this tool is useful?

What are the constrains or difficulties on implementing this model for reliable analysis?

VPD in (kPa)

vpd =ei-ea_v round(vpd, digits=4 )